# 渗流模型动画程序

# 程序分为五大部分
	# 一、库调用模块。
	# 二、初始参数设定模块
	# 三、团簇计算模块
	# 四、主函数更新模块
	# 五、动画生成模块


# 一、库调用模块。
# 调用库函数
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
import random
import sys


# 设置递归深度为100000
# 因为本程序多次使用for循环和嵌套递归，默认情况下递归深度达到1000会出现异常
sys.setrecursionlimit(100000) #这里设置为十万


# 创建画布
fig = plt.figure()


# 二、初始参数设定模块
# 初始值设定
# 确定总体渗流模型的规模
# 最开始的渗流模型是全 1 矩阵
# Bao是一比一对应渗流模型坐标的坐标矩阵
# Bao矩阵作用是为了防止重复选点
Bao = []
boardsize = int(input('请输入边长:'))
my_board = np.ones((boardsize , boardsize ))
for i in range(boardsize):
	for j in range(boardsize):
		Bao.append([i, j])


# 生成全0 N矩阵
# 该矩阵是为了统计在几行几列的位置上存在最大、次大团簇
N = np.zeros((boardsize , boardsize))


# 设置团簇的生成概率P
P = []
for i in range(boardsize*boardsize):
	P.append((i+4)/(boardsize*boardsize)) # 因为程序的遍历原因，初始值点会从3起步，加上i从0开始。所以i+4


# 设置最大团簇计算矩阵（N）索引的初值
l = 0  # l是上面生成的矩阵N行数的索引控制值
k = -1  # k是上面生成的矩阵N列数的索引控制值

# 三、团簇计算模块
# 方向控制遍历函数
# 该部分一共定义了4*3=12个方向控制遍历函数（up、down、left、right）,每一个方向又可分为序号1、2、3三个不同功能的函数
# 1：负责统计当前矩阵有多少0值，即已经生成了多少个点
# 2：判断哪个团簇为最大团簇
# 3：判断哪个团簇为次大团簇
def up(i, j):  # 负责控制团簇向上方向遍历
	global l
	global k
	if i == -boardsize-1:	# 因为矩阵的遍历范围只有（-boardsize,boardsize-1）为了防止遍历时超出索引范围，需要重置索引值
		i = boardsize-1
	if new_arr[i][j] == 0: # 当遍历到的值为0时，说明该点为团簇
		N[l][k] = N[l][k] + 1	# 每发现一个团簇便记录一次团簇大小
		new_arr[i][j] = 2	#每当记录完一个团簇便重置该点的值为2，以防之后的重复计算
		e = i
		f = j - 1
		left(e, f)
		g = i
		h = j + 1
		right(g, h)
		i = i - 1
		up(i, j)
	else:
		N[l][k] = N[l][k]


def up2(i, j):  # 负责控制团簇向上方向遍历
	if i == -boardsize-1:
		i = boardsize-1
	if new_arr[i][j] == 2:
		new_arr[i][j] = 3
		e = i
		f = j - 1
		left2(e, f)
		g = i
		h = j + 1
		right2(g, h)
		i = i - 1
		up2(i, j)


def up3(i, j):  # 负责控制团簇向上方向遍历
	if i == -boardsize-1:
		i = boardsize-1
	if new_arr[i][j] == 2:
		new_arr[i][j] = 4
		e = i
		f = j - 1
		left3(e, f)
		g = i
		h = j + 1
		right3(g, h)
		i = i - 1
		up3(i, j)


def down(i, j):  # 负责控制团簇向下方向遍历
	global k
	global l
	if i == boardsize:
		i = -boardsize
	if new_arr[i][j] == 0:
		N[l][k] = N[l][k] + 1
		new_arr[i][j] = 2
		e = i
		f = j - 1
		left(e, f)
		g = i
		h = j + 1
		right(g, h)
		i = i + 1
		down(i, j)
	else:
		N[l][k] = N[l][k]


def down2(i, j):  # 负责控制团簇向下方向遍历
	if i == boardsize:
		i = -boardsize
	if new_arr[i][j] == 2:
		new_arr[i][j] = 3
		e = i
		f = j - 1
		left2(e, f)
		g = i
		h = j + 1
		right2(g, h)
		i = i + 1
		down2(i, j)


def down3(i, j):  # 负责控制团簇向下方向遍历
	if i == boardsize:
		i = -boardsize
	if new_arr[i][j] == 2:
		new_arr[i][j] = 4
		e = i
		f = j - 1
		left3(e, f)
		g = i
		h = j + 1
		right3(g, h)
		i = i + 1
		down3(i, j)


def left(i, j):  # 负责控制团簇向左方向遍历
	global k
	global l
	if j == -boardsize-1:
		j = boardsize-1
	if new_arr[i][j] == 0:
		N[l][k] = N[l][k] + 1
		new_arr[i][j] = 2
		a = i - 1
		b = j
		up(a, b)
		c = i + 1
		d = j
		down(c, d)
		j = j - 1
		left(i, j)
	else:
		N[l][k] = N[l][k]


def left2(i, j):  # 负责控制团簇向左方向遍历
	if j == -boardsize-1:
		j = boardsize-1
	if new_arr[i][j] == 2:
		new_arr[i][j] = 3
		a = i - 1
		b = j
		up2(a, b)
		c = i + 1
		d = j
		down2(c, d)
		j = j - 1
		left2(i, j)


def left3(i, j):  # 负责控制团簇向左方向遍历
	if j == -boardsize-1:
		j = boardsize-1
	if new_arr[i][j] == 2:
		new_arr[i][j] = 4
		a = i - 1
		b = j
		up3(a, b)
		c = i + 1
		d = j
		down3(c, d)
		j = j - 1
		left3(i, j)


def right(i, j):  # 负责控制团簇向右方向遍历
	global k
	global l
	if j == boardsize:
		j = -boardsize
	if new_arr[i][j] == 0:
		N[l][k] = N[l][k] + 1
		new_arr[i][j] = 2
		a = i - 1
		b = j
		up(a,b)
		c = i + 1
		d = j
		down(c,d)
		j = j + 1
		right(i,j)
	else:
		N[l][k] = N[l][k]


def right2(i, j):  # 负责控制团簇向右方向遍历
	if j == boardsize:
		j = -boardsize
	if new_arr[i][j] == 2:
		new_arr[i][j] = 3
		a = i - 1
		b = j
		up2(a, b)
		c = i + 1
		d = j
		down2(c, d)
		j = j + 1
		right2(i, j)


def right3(i, j):  # 负责控制团簇向右方向遍历
	if j == boardsize:
		j = -boardsize
	if new_arr[i][j] == 2:
		new_arr[i][j] = 4
		a = i - 1
		b = j
		up3(a, b)
		c = i + 1
		d = j
		down3(c, d)
		j = j + 1
		right3(i, j)


# 团簇计算模块
# 最大团簇计算模块
# 通过找到最大团簇的坐标max，以此为起点开始遍历该团簇，将该团簇的值赋予3
def change_color(max):
	i = max[0][0]
	j = max[1][0]
	new_arr[i][j] = 3
	a = i - 1
	b = j
	up2(a, b)
	c = i + 1
	d = j
	down2(c, d)
	e = i
	f = j - 1
	left2(e, f)
	g = i
	h = j + 1
	right2(g, h)


# 次大团簇计算模块
# 通过找到次大团簇的坐标max2，以此为起点开始遍历该团簇，将该团簇的值赋予4
def change_color2(max2):
	i = max2[0][0]
	j = max2[1][0]
	new_arr[i][j] = 4
	a = i - 1
	b = j
	up3(a, b)
	c = i + 1
	d = j
	down3(c, d)
	e = i
	f = j - 1
	left3(e, f)
	g = i
	h = j + 1
	right3(g, h)


def Sagiri(my_board):
	global N
	global k
	global l
	global new_arr
	global max,max2

	# 随机选点赋值
	bao = random.choice(Bao)
	my_board[bao[0]][bao[1]] = 0
	Bao.remove(bao)	#该函数能将已经选择的坐标移出Bao列表，防止重复选点

	for i in range(0, boardsize, 1):
		for j in range(0, boardsize, 1):
			if my_board[i][j] == 3 or my_board[i][j] == 2 or my_board[i][j] == 4 or my_board[i][j] == 6:
				my_board[i][j] = 0


	new_arr = my_board


	for i in range(0, boardsize, 1):
		for j in range(0, boardsize, 1):
			l = i
			k = j
			if new_arr[i][j] == 0:
				N[l][k] = N[l][k] + 1
				new_arr[i][j] = 2
				a = i - 1
				b = j
				up(a, b)
				c = i + 1
				d = j
				down(c, d)
				e = i
				f = j - 1
				left(e, f)
				g = i
				h = j + 1
				right(g, h)

			if new_arr[i][j] == 1:  # 当矩阵中值为1时，不计入团簇，并结束该次循坏
				continue
	max = np.where(N == np.max(N))

	change_color(max)

	# 暂存最大团簇的坐标值，并暂时消除最大团簇的影响，方便计算次大团簇
	for i in range(0, boardsize, 1):
		for j in range(0, boardsize, 1):
			if my_board[i][j] == 2:
				my_board[i][j] = 0

	N = np.zeros((boardsize , boardsize ))

	for i in range(0, boardsize, 1):
		for j in range(0, boardsize, 1):
			l = i
			k = j
			if new_arr[i][j] == 0:
				N[l][k] = N[l][k] + 1
				new_arr[i][j] = 2
				a = i - 1
				b = j
				up(a, b)
				c = i + 1
				d = j
				down(c, d)
				e = i
				f = j - 1
				left(e, f)
				g = i
				h = j + 1
				right(g, h)

			if new_arr[i][j] == 1:  # 当矩阵中值为1时，不计入团簇，并结束该次循坏
				continue

	max2 = np.where(N == np.max(N))#由于第一次执行该循环时，N为全0矩阵。会导致默认次大团簇为[0,0]，所以会导致[0,0]在每次程序运行时，第一次循环予其固定赋予值4。
									# 但由于修改此处会导致BUG，且影响不大，故不做修改。仅在此处予以说明。
	change_color2(max2)


# 四、主函数更新模块
def update_board(my_board):
	global N
	global k
	global l
	global new_arr
	global max

	Sagiri(my_board)
	N = np.zeros((boardsize , boardsize ))
	new_board = new_arr
	return new_board


# 五、动画生成模块


update_board(my_board)# 具体运作原理不明，但没有该行无法画图
						# 猜测原理是在动画开始前先给出一个初始矩阵

# imshow()函数可以接收一个numpy矩阵并返回一张图片
im = plt.imshow(my_board, animated=True)


# 负责控制动画
# animation是作为FuncAnimation() 调用的辅助函数。
# animate() 函数接受一帧画面作为输入充当计数器。
# 这个画面计数器就是  FuncAnimation() 和 animate() 函数沟通的桥梁——在每一个时间点（也就是每一帧），它都会调用一次  animate()。
# 然后 animate() 会逐次使用辅助函数  update_board() 来对网格进行迭代。最后，set_data() 函数将图片更新为迭代后的网格。
def animate(num):
	im.set_array(update_board(my_board))  # 用来更新矩阵中的值，以此实现更新矩阵中色块的目的
	plt.title('P=%.4f'%P[num])	# 负责更新概率P
	return im


# 调用刚刚定义好的动画函数
# frames：代表画面的帧数
# fig：是我们在前面创建的用来装载我们的动画的图形变量；
# animate：是FuncAnimation()用画面计数器进行沟通的函数（自动传入，不需要特别声明）
# frames：表示我们希望动画持续多少帧，在这里我们想要动画的长度为多少帧；
# interval：表示每一帧之间间隔的毫秒数。我们想要每帧之间间隔多少毫秒。
anim = animation.FuncAnimation(fig, func=animate, frames=50000, interval=0.0001, blit=False)


# 该函数可以生成gif图
anim.save('Baos（100*100）.gif')


# 绘制图像
plt.show()